#!/bin/tcsh -ef
####
#
#############################################################################
#                                                                           #
# Title: Merge to Make Reference                                            #
#                                                                           #
# (C) 2dx.org, GNU Plublic License.                                         #
#                                                                           #
# Created..........: 01/03/2007                                             #
# Last Modification: 01/03/2007                                             #
# Author...........: 2dx.org                                                #
#                                                                           #
#############################################################################
#
# SORTORDER: 15
#
# MANUAL: This script allows merging the selected image data into one merge dataset.
#
# MANUAL: Select only the one best non-tilted image and run this script with that one image.
#
# MANUAL: Then align the other images to it (using Refine Once), starting with non-tilted images, and slowly add more and more tilted images. After each call of "Refine Once" you should run "Merge Once" again, to merge those aligned images into a new reference. 
#
# MANUAL: NOTE: Select ONLY the non-tilted images and run this script!
#
#
# PUBLICATION: 3D reconstruction of two-dimensional crystals: <A HREF="http://www.ncbi.nlm.nih.gov/pubmed/26093179">Arch Biochem Biophys 581, 68-77 (2015)</A>
# PUBLICATION: 3D Reconstruction from 2D crystal image and diffraction data: <A HREF="http://dx.doi.org/10.1016/S0076-6879(10)82004-X">Methods Enzymol. 482, Chapter 4, 101-129 (2010)</A>
# PUBLICATION: 2dx - Automated 3D structure reconstruction from 2D crystal data: <A HREF="http://journals.cambridge.org/action/displayAbstract?aid=1943200">Microscopy and Microanalysis 14(Suppl. 2), 1290-1291 (2008)</A>
# PUBLICATION: 2dx_merge - Data management and merging for 2D crystal images: <A HREF="http://dx.doi.org/10.1016/j.jsb.2007.09.011">J. Struct. Biol. 160(3), 375-384 (2007)</A>
#
# DISPLAY: MergeResolution
# DISPLAY: RESMIN
# DISPLAY: RESMAX
# DISPLAY: merge_res_limit
# DISPLAY: merge_data_type
# DISPLAY: tempkeep
# DISPLAY: realcell
# DISPLAY: realang
# DISPLAY: ALAT
# DISPLAY: SYM
# DISPLAY: avrgamphsNUMBER
# DISPLAY: avrgamphsRESOL
# DISPLAY: ILIST
# DISPLAY: MergeIQMAX
# DISPLAY: Merge_Reference_IQMAX
# DISPLAY: MergeHKMAX
# DISPLAY: Merge_Reference_HKMAX
# DISPLAY: RFACAMP
# DISPLAY: MergeAK
# DISPLAY: MergeIWF
# DISPLAY: AMP_Scale_Factor
# DISPLAY: MergeIWP
# DISPLAY: max_amp_correction
# DISPLAY: plotres_rings
# DISPLAY: resolutionplot_RESMAX
# DISPLAY: resolutionplot_bins
# DISPLAY: Thread_Number
# DISPLAY: RESMIN
# DISPLAY: RESMAX
#
#$end_local_vars
#
set bin_2dx = ""
set proc_2dx = ""
#
set ALAT = ""
set AMP_Scale_Factor = ""
set avrgamphsNUMBER = ""
set avrgamphsRESOL = ""
set ILIST = ""
set lattice = ""
set max_amp_correction = ""
set merge_data_type = ""
set merge_ref_num = ""
set Merge_Reference_HKMAX = ""
set Merge_Reference_IQMAX = ""
set merge_res_limit = ""
set MergeAK = ""
set MergeHKMAX = ""
set MergeIQMAX = ""
set MergeIWF = ""
set MergeIWP = ""
set MergeResolution = ""
set num_amplitudes_observed = ""
set num_phases_observed = ""
set num_reflections_fitted = ""
set num_reflections_FOM50 = ""
set overall_phase_residual = ""
set overall_R_factor = ""
set overall_weighted_phase_residual = ""
set overall_weighted_R_factor = ""
set plotres_rings = ""
set realang = ""
set realcell = ""
set Reflections_Unique = ""
set RESMAX = ""
set RESMIN = ""
set resolutionplot_bins = ""
set resolutionplot_RESMAX = ""
set RFACAMP = ""
set sample_pixel = ""
set SYM = ""
set tempkeep = ""
set Thread_Number = ""
set zstarwin = ""
set ctfrev = ""
#
#$end_vars
#
# 
set scriptname = 2dx_merge2D
#
set merge_modus="2D"
#
\rm -f LOGS/${scriptname}.results
#
set ccp4_setup = 'y'
source ${proc_2dx}/initialize
#
# The following is to make sure that for the next "Import Images", the default is correctly initialized.
set initialization_reset = "y"
set initialization_executable = "y"
echo "set initialization_reset = ${initialization_reset}" >> LOGS/${scriptname}.results
echo "set initialization_executable = ${initialization_executable}" >> LOGS/${scriptname}.results
#
#
#----------------------------------------------------------------------
${proc_2dx}/linblock "Verifying some parameters"
#----------------------------------------------------------------------
#
if ( `echo ${RESMIN} ${RESMAX} | awk '{ if ( $1 < $2 ) { s = 1 } else { s = 0 }} END { print s }'` == 1 ) then
    set oldval = ${RESMIN}
    set RESMIN = ${RESMAX}
    set RESMAX = ${oldval}
    ${proc_2dx}/linblock "ERROR: exchanging RESMIN and RESMAX, to RESMIN=${RESMIN}, and RESMAX=${RESMAX}"
    echo "set RESMIN = ${RESMIN}" >> LOGS/${scriptname}.results
    echo "set RESMAX = ${RESMAX}" >> LOGS/${scriptname}.results
endif
#
if ( ${merge_data_type}x == 'x' ) then
    set merge_data_type = "0"
    echo "set merge_data_type = ${merge_data_type}" >> LOGS/${scriptname}.results
    ${proc_2dx}/linblock ":WARNING: merge_data_type corrected to ${merge_data_type}"
endif
#
if ( ${MergeHKMAX}x == 'x' ) then
    set MergeHKMAX = '20'
    ${proc_2dx}/linblock "ERROR: correcting MergeHKMAX to ${MergeHKMAX}"
    echo "set MergeHKMAX = ${MergeHKMAX}" >> LOGS/${scriptname}.results
endif
#
#
if ( ${AMP_Scale_Factor}x == 'x' ) then
    set AMP_Scale_Factor = 0.001
    echo "set AMP_Scale_Factor = ${AMP_Scale_Factor}" >> LOGS/${scriptname}.results
    ${proc_2dx}/linblock "WARNING: AMP_Scale_Factor corrected to ${AMP_Scale_Factor}"
endif
# 
if ( ${MergeAK}x == 'x' ) then
    set MergeAK = "10"
    echo "set MergeAK = ${MergeAK}" >> LOGS/${scriptname}.results
    ${proc_2dx}/linblock "WARNING: MergeAK corrected to ${MergeAK}"
endif
#
if ( ${merge_res_limit}x == 'x' ) then
    set merge_res_limit = "0"
    echo "set merge_res_limit = ${merge_res_limit}" >> LOGS/${scriptname}.results
    ${proc_2dx}/linblock "WARNING: merge_res_limit corrected to ${merge_res_limit}"
endif
#
if ( ${MergeIWF}x == 'x' ) then
    set MergeIWF = "0"
    echo "set MergeIWF = ${MergeIWF}" >> LOGS/${scriptname}.results
    ${proc_2dx}/linblock "WARNING: MergeIWF corrected to ${MergeIWF}"
endif
#
if ( ${MergeIWP}x == 'x' ) then
    set MergeIWP = "1"
    echo "set MergeIWP = ${MergeIWP}" >> LOGS/${scriptname}.results
    ${proc_2dx}/linblock "WARNING: MergeIWP corrected to ${MergeIWP}"
endif
#
if ( ${ILIST}x == 'x' ) then
    set ILIST = "n"
    echo "set ILIST = ${ILIST}" >> LOGS/${scriptname}.results
    ${proc_2dx}/linblock "WARNING: ILIST corrected to ${ILIST}"
endif
#
if ( ${ILIST} == "n" ) then
    set ILIST_VAL = 0
else
    set ILIST_VAL = 1
endif
#
if(${MergeIWF} == "0")then
    set MergeIWF_VAL = -1
endif
if(${MergeIWF} == "1")then
    set MergeIWF_VAL = 0
endif
if(${MergeIWF} == "2")then
    set MergeIWF_VAL = 1
endif
if(${MergeIWP} == "0")then
    set MergeIWP_VAL = -2
endif
if(${MergeIWP} == "1")then
    set MergeIWP_VAL = -1
endif
if(${MergeIWP} == "2")then
    set MergeIWP_VAL = 0
endif
if(${MergeIWP} == "3")then
    set MergeIWP_VAL = 1
endif
#
echo ILIST_VAL = ${ILIST_VAL}
echo MergeAK = ${MergeAK}
echo MergeIWF_VAL = ${MergeIWF_VAL}
echo MergeIWP_VAL = ${MergeIWP_VAL}
#
# This memorizes the current merge directory under the variable "olddir":
set olddir = $PWD
#
echo "The current working directory is" ${olddir}
#
#----------------------------------------------------------------------
${proc_2dx}/linblock "Sourcing ${proc_2dx}/2dx_check_links.com ${olddir}"
#----------------------------------------------------------------------
source ${proc_2dx}/2dx_check_links.com ${olddir}
#
#
#This translates the list of directories to work on into one single long line:
cat 2dx_merge_dirfile.dat | tr "\n" " " > SCRATCH/2dx_merge_dirfile_oneline.dat
set dirlist = "`cat SCRATCH/2dx_merge_dirfile_oneline.dat`"
#
#
#---------------------------------------------------------------------------------------------------
${proc_2dx}/linblock "Sourcing ${proc_2dx}/2dx_select_APH.com  ${olddir} ${merge_data_type}"
#----------------------------------------------------------------------------------------------------
source ${proc_2dx}/2dx_select_APH.com  ${olddir} ${merge_data_type}
#
#
foreach dirfile ( ${dirlist} ) 
    cd ..
    cd ${dirfile}
    set imagename_local = `cat 2dx_image.cfg | grep 'set imagename =' | cut -d\" -f2`
    cd APH
    if ( ! -e image_ctfcor_fou_unbent_ctf.aph ) then
        if ( -e ${imagename_local}.cor.aph ) then
            echo "::In ${dirfile}: Correcting ${imagename_local}.cor.aph to image_ctfcor_fou_unbent_ctf.aph" 
            \mv -f ${imagename_local}.cor.aph image_ctfcor_fou_unbent_ctf.aph 
        else
            echo "::WARNING: ${dirfile}/APH/image_ctfcor_fou_unbent_ctf.aph does not exist. "
        endif
    else
        # echo "Working on ${dirfile}/APH/image_ctfcor_fou_unbent_ctf.aph exists. "
    endif
    cd ${olddir}
end
#
echo "<<@progress: 5>>"
#
#----------------------------------------------------------------------
${proc_2dx}/linblock "sourcing sym2spsgrp_sub.com"
#----------------------------------------------------------------------
#
source ${proc_2dx}/2dx_sym2spcgrp_sub.com
#
echo SYM = ${SYM}
echo spcgrp = ${spcgrp}
echo CCP4_SYM = ${CCP4_SYM}
#
#----------------------------------------------------------------------
${proc_2dx}/lin "2dx_merge_makedirs - to create all required subdirectories"
#----------------------------------------------------------------------
#
source ${proc_2dx}/2dx_merge_makedirs
#
echo "<<@progress: 10>"
#
#----------------------------------------------------------------------
${proc_2dx}/linblock "Sourcing 2dx_origtilt_merge.com "
#----------------------------------------------------------------------
#
set dirfile = "2dx_merge_dirfile.dat"
set select_mode = "all"
#
source ${proc_2dx}/2dx_origtilt_merge.com
#
source ${proc_2dx}/2dx_origtilt_aph_to_mtz.com
#----------------------------------------------------------------------
${proc_2dx}/linblock "${scriptname} normal end."
#----------------------------------------------------------------------
#
echo "<<@progress: 100>>"
#
